Dans un premier temps, nous utilisons la fonction « list.files » pour récupérer le nom de tous les fichiers débutant par « bnvd », « bnvd » étant le début commun à toutes les bases de données concernant l’achats de produits phytosanitaires.
Ensuite, nous retirons de cette base de données toutes les lignes où l’identifiant est « 0000 », puis les valeurs manquantes sont transformées en « NA » (ce qui n’était pas le cas avant).
Après, nous enlevons les lignes pour lesquelles la quantité de substance n’est pas renseignée puis, enfin, nous regroupons certains facteurs sous le nom « substance ».
Par la suite, un premier « group_by » est effectué. Pour un même code postal, la même année et la même classification, nous additionnons la quantité de substances achetées.
Le deuxième « group_by » reprend le dataframe obtenu lors du premier « group_by » et permet, cette fois, d’obtenir une seule ligne pour un code postal.
fichiers = list.files(pattern="^bnvd.*\\.csv$")
bddsubstance = do.call(rbind, lapply(fichiers, function(x) read.csv(x, sep=";", dec=".", stringsAsFactors = FALSE,header=TRUE, colClasses=c("code_postal_acheteur"="character"),encoding = "UTF-8" )))
bddsubstance <- bddsubstance[-which(bddsubstance$code_postal_acheteur=="00000"),]
bddsubstance$quantite_substance[bddsubstance$quantite_substance == ""] <- NA
bddsubstance<-bddsubstance[!is.na(bddsubstance$quantite_substance),]
bddsubstance$quantite_substance<-as.numeric(bddsubstance$quantite_substance)
## Warning: NAs introduits lors de la conversion automatique
bddsubstance$classification<-as.factor(bddsubstance$classification)
bddsubstance<-bddsubstance[c(1,2,6,7)]
levels(bddsubstance$classification)[c(4,5)]<-"Toxique"
bddsubstance = bddsubstance %>% group_by(code_postal_acheteur,annee, classification) %>%
summarise(quantite_totale = sum(quantite_substance), .groups = 'drop')
bddlongsubstance<-bddsubstance %>%
pivot_wider(
id_cols = code_postal_acheteur,
names_from = c(annee, classification),
values_from = quantite_totale,
names_glue = "{classification}_quantite-totale-{annee}"
)
Voici le résultat du premier « group_by ».
Nous obtenons une quantité totale unique par code postal, par année et par classification.
## # A tibble: 85,752 x 4
## code_postal_acheteur annee classification quantite_totale
## <chr> <int> <fct> <dbl>
## 1 01000 2015 Autre 44.3
## 2 01000 2015 N minéral 30.2
## 3 01000 2015 N Organique 690.
## 4 01000 2015 Toxique 213.
## 5 01000 2016 Autre 386.
## 6 01000 2016 N minéral 479.
## 7 01000 2016 N Organique 2671.
## 8 01000 2016 Toxique 274.
## 9 01000 2017 Autre 1460.
## 10 01000 2017 N minéral 275.
## # ... with 85,742 more rows
Le deuxième « group_by » permet de n’avoir qu’une seule ligne par code postal, le transformant en identifiant.
## # A tibble: 6,126 x 17
## code_postal_acheteur `Autre_quantite-tota~` `N minéral_qua~` `N Organique_q~`
## <chr> <dbl> <dbl> <dbl>
## 1 01000 44.3 30.2 690.
## 2 01090 164. 7.60 2715.
## 3 01100 0.450 2 145.
## 4 01110 3.85 4.3 132.
## 5 01120 328. 14.3 8800.
## 6 01130 0.00745 NA 24.2
## 7 01140 392. 21.4 7377.
## 8 01150 314. 69.4 7123.
## 9 01160 115. 20.1 2391.
## 10 01170 396. 4.76 1119.
## # ... with 6,116 more rows, and 13 more variables:
## # `Toxique_quantite-totale-2015` <dbl>, `Autre_quantite-totale-2016` <dbl>,
## # `N minéral_quantite-totale-2016` <dbl>,
## # `N Organique_quantite-totale-2016` <dbl>,
## # `Toxique_quantite-totale-2016` <dbl>, `Autre_quantite-totale-2017` <dbl>,
## # `N minéral_quantite-totale-2017` <dbl>,
## # `N Organique_quantite-totale-2017` <dbl>, ...
De la même manière que pour la base de données sur les quantités achetées, nous utilisons la fonction « list.files » pour détecter tous les fichiers dont le nom commence par « moy ».
Le dataframe contient les valeurs qui nous intéressent, cependant l’identifiant est « CD_station ». Or, nous souhaitons le code postal.
Pour cela, nous importons « stat », qui nous donne, entre autres, la station (avec son identifiant « CD_station »), son code postal, son code Insee.
Le dataframe « cp », lui, associe à un code Insee le code postal correspondant. Par le biais de deux « merge », nous parvenons à obtenir un dataframe avec les données qui nous intéressent, et le code postal en identifiant.
Enfin, un « group_by » est effectué de manière à avoir le code postal comme identifiant.
fichiers = list.files(pattern="^moy.*\\.csv$")
pe<- do.call(rbind, lapply(fichiers, function(x) read.csv(x, sep=";", dec=",", stringsAsFactors = FALSE,header=TRUE,encoding = "UTF-8" )))
cp<-read.table(file="codepostal.csv",header=TRUE, sep=";", dec=',',colClasses=c("Code_postal"="character"),encoding = "UTF-8")
stat<-read.csv2(file="stations.csv")
stat<-stat[,c(1,2)]
cp<-cp[,c(1,3)]
fus<-merge(pe,stat, by="CD_STATION")
names(cp)<-c("NUM_COM","COD_POS")
bddeau<-merge(fus,cp,by="NUM_COM")
bddeau<-bddeau[,-c(1,2)]
bddlongeau = bddeau %>% group_by(COD_POS,ANNEE) %>%
summarise(nbpreltot=sum(NBPREL),MOYPTOTAL=mean(MOYPTOT), MAXPTOTAL=max(MAXPTOT), MINMOLRECHTOTAL=min(MINMOLRECH), MINMOLQTOTAL=min(MINMOLQ),MAXMOLQTOTAL=max(MAQMOLQ), .groups = 'drop')
rm(list=c("cp","fus","pe","stat","fichiers"))
Voici les deux dataframes obtenus :
## ANNEE NBPREL MOYPTOT MAXPTOT MINMOLRECH MAXMOLRECH MINMOLQ MAQMOLQ COD_POS
## 1 2007 6 0.05333333 0.09 376 377 1 3 01500
## 2 2010 5 0.32900000 0.56 247 409 2 3 01360
## 3 2007 7 0.10028571 0.18 298 377 1 3 01360
## 4 2011 4 0.25000000 0.36 409 409 2 3 01360
## 5 2012 5 0.21720000 0.33 247 409 2 2 01360
## 6 2009 5 0.26260000 0.43 299 383 2 3 01360
## # A tibble: 6 x 8
## COD_POS ANNEE nbpreltot MOYPTOTAL MAXPTOTAL MINMOLRECHTOTAL MINMOLQTOTAL
## <chr> <int> <int> <dbl> <dbl> <int> <int>
## 1 01000 2009 8 0 0 2 0
## 2 01000 2010 4 0.04 0.06 370 1
## 3 01000 2011 8 0.03 0.04 335 0
## 4 01000 2012 8 0.025 0.04 368 1
## 5 01120 2007 2 0 0 1 0
## 6 01120 2008 4 0 0 382 0
## # ... with 1 more variable: MAXMOLQTOTAL <int>
On transforme le code postal en facteur, puis nous remplaçons « code_postal_acheteur », le nom de la colonne contenant les codes postaux, par « ID ».
Enfin, nous importons le fichier shapefile codes_postaux_région.shp avant d’encoder le fichier en UTF-8.
bddlongsubstance$code_postal_acheteur= as.factor(bddlongsubstance$code_postal_acheteur)
names(bddlongsubstance)[1]= "ID"
codes_postaux= st_read(dsn = "codes_postaux_V5/codes_postaux_region.shp",
layer = "codes_postaux_region",
quiet = TRUE) %>%
select(ID, LIB, DEP)
codes_postaux$LIB= stri_encode(codes_postaux$LIB, from= "ISO-8859-1", to= "utf8")
Nous fusionnons ensuite le fichier shapefile avec la base de données.
fusion= codes_postaux %>%
left_join(bddlongsubstance[1:17], by= "ID") %>%
st_transform(2154)
Enfin, nous créons la carte pour les valeurs « autre » de l’année 2015.
Ce processus est le même pour toutes les autres catégories, il suffit juste de changer le paramètre « col » dans « tm_fill » et « quantité » dans « popup.vars ».
tmap_mode(mode= "view")
## tmap mode set to interactive viewing
autre_2015=
tm_basemap("CartoDB.Voyager")+ tm_shape(shp= fusion)+
tm_fill(col= "Autre_quantite-totale-2015", palette= "YlOrBr", id= "ID",
textNA = "Valeur manquante", style = "quantile", n= 6, title= "Quantité totale : <br> Catégorie AUTRE",
popup.vars = c("Ville"= "LIB", "Quantité"= "Autre_quantite-totale-2015"))+
tm_borders("black", lwd= 0.3, alpha= 0.6)+
tm_layout(title = "Quantité de substances phytopharmaceutiques achetées en 2015 (en kilogrammes)", title.position = c("center", "bottom"), legend.bg.color = "white", legend.bg.alpha = 0.4)+
tm_scale_bar(position = c("left", "bottom"))+
tm_view(view.legend.position = c("right", "bottom"))
Nous affichons ensuite la carte créée pour la catégorie « autre » en 2015.